#include "VTECConverter.h"
#include <cmath>
#include <iostream>

// 常量定义
constexpr double PI = 3.14159265358979323846;
constexpr double DEG2RAD = PI / 180.0;
constexpr double EARTH_RADIUS = 6371.0; // 地球平均半径(km)

VTECConverter::VTECConverter(double ionosphereHeight) 
    : ionosphereHeight_(ionosphereHeight) {
    // 初始化
}

VTECConverter::~VTECConverter() {
    // 清理资源
}

void VTECConverter::setSTECObservations(const std::vector<STECObservation>& stecObs) {
    stecObservations_ = stecObs;
}

void VTECConverter::convertToVTEC() {
    vtecObservations_.clear();
    
    for (const auto& stecObs : stecObservations_) {
        // 计算映射函数
        double mappingFunc = calculateMappingFunction(stecObs.elevation);
        
        // 从STEC计算VTEC
        double vtec = stecObs.stec * mappingFunc;
        
        // 创建VTEC观测对象
        VTECObservation vtecObs;
        vtecObs.time = stecObs.time;
        vtecObs.satellite = stecObs.satellite;
        vtecObs.station = stecObs.station;
        vtecObs.latitude = stecObs.latitude;
        vtecObs.longitude = stecObs.longitude;
        vtecObs.elevation = stecObs.elevation;
        vtecObs.stec = stecObs.stec;
        vtecObs.vtec = vtec;
        vtecObs.mappingFunction = mappingFunc;
        
        vtecObservations_.push_back(vtecObs);
    }
    
    std::cout << "转换得到 " << vtecObservations_.size() << " 个VTEC观测值" << std::endl;
}

std::vector<VTECObservation> VTECConverter::getVTECObservations() const {
    return vtecObservations_;
}

double VTECConverter::calculateMappingFunction(double elevation) {
    // 转换为弧度
    double elev = elevation * DEG2RAD;
    
    // 计算映射函数: 1/cos(z')
    // 其中z'是穿刺点的天顶角
    double z = PI/2 - elev;
    double zPrime = std::asin((EARTH_RADIUS / (EARTH_RADIUS + ionosphereHeight_)) * std::sin(z));
    
    return 1.0 / std::cos(zPrime);
} 